Analysis of stimulation experiments.

library(Seurat)
## Attaching SeuratObject
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(ggplot2)

Load integrated object. This object has been integrated with k.anchor = 20, because there was lots of batch effect between the 0hr tissue and the other two experiments (which merged quite well).

load("~/Dropbox/Zoe/scf_version/analysis/healthy_sc/seurat_objects/no_dropletQC/3391/integrated_3391_stimulation_kanchor20_scClustViz.RData")

Set identities:

Idents(scSeurat) <- "integrated_snn_res.0.4"

Check out the UMAP of clusters, original identities, and refined SCINA labels.

DimPlot(scSeurat, label = TRUE)

DimPlot(scSeurat, group.by = "tissue_type", label = FALSE)

DimPlot(scSeurat, group.by = "scina_labels_refined", label = TRUE) & NoLegend()

Look at tissue-type per cluster on a proportional level:

pt <- table(Idents(scSeurat), scSeurat$tissue_type)
pt <- as.data.frame(pt)
ggplot(pt, aes(x = Var1, y = Freq, fill = Var2)) +
  theme_bw(base_size = 15) +
  geom_col(position = "fill", width = 0.5) +
  xlab("Sample") +
  ylab("Proportion") +
  theme(legend.title = element_blank())

Look for T-cell signatures.

FeaturePlot(scSeurat, features = c("CD3E", "NKG7"))

FeaturePlot(scSeurat, features = c("LEF1", "XCL1;XCL2"))

FeaturePlot(scSeurat, features = c("TOX", "CD8A"))

Look for macrophage signatures:

FeaturePlot(scSeurat, features = c("LGALS1", "TYROBP"))

FeaturePlot(scSeurat, features = c("MARCO", "C1QC"))

FeaturePlot(scSeurat, features = c("LYZ-1", "S100A8"))

FeaturePlot(scSeurat, features = c("DEFA6;DEFA4;DEFA5", "CTSS"))

FeaturePlot(scSeurat, features = c("FLT3", "S100A12"))

Cluster 11 is the only cluster that is obviously dominated by the 0hr and unstimulated experiments. What cells are these? Looks like they are hepatocytes. A lot of the markers are periportal, but all hepatocytes are in this cluster.

# Isolate DE genes
clust11 <- sCVdata_list$res.0.4@DEvsRest$`11`[order(sCVdata_list$res.0.4@DEvsRest$`11`$Wstat,
                                                    decreasing = TRUE), ]
rownames(clust11)[1:50]
##  [1] "PAH"                                     
##  [2] "LOC114082496"                            
##  [3] "GYS2"                                    
##  [4] "SLC2A2"                                  
##  [5] "SLCO1B3;SLCO1B3-SLCO1B7;SLCO1B7;SLCO1B1" 
##  [6] "PCK1"                                    
##  [7] "CPS1"                                    
##  [8] "AGMO"                                    
##  [9] "MTSS1"                                   
## [10] "FBP1"                                    
## [11] "ALDOB"                                   
## [12] "MAT1A"                                   
## [13] "mikado.WCK01-AAH20201022-F8-rsc00143G668"
## [14] "SHMT1"                                   
## [15] "CYP3A4;CYP3A7-CYP3A51P;CYP3A7-3"         
## [16] "KLKB1"                                   
## [17] "TDO2"                                    
## [18] "SLC38A4"                                 
## [19] "SORBS2"                                  
## [20] "SOX5"                                    
## [21] "CYP3A4;CYP3A7-CYP3A51P;CYP3A7-2"         
## [22] "ARG1"                                    
## [23] "SLC25A13"                                
## [24] "CYP3A4;CYP3A7-CYP3A51P;CYP3A7"           
## [25] "Slc22a19-1"                              
## [26] "TAT"                                     
## [27] "LOC107139912"                            
## [28] "KIRREL3"                                 
## [29] "ACACB"                                   
## [30] "C5-1"                                    
## [31] "SUGCT"                                   
## [32] "ZBTB16"                                  
## [33] "ABCB11"                                  
## [34] "CFH-1"                                   
## [35] "GAS2"                                    
## [36] "ACSL5"                                   
## [37] "STARD13"                                 
## [38] "IGFBP1"                                  
## [39] "HSD11B1"                                 
## [40] "AQP9"                                    
## [41] "GRHPR"                                   
## [42] "FMO5-1"                                  
## [43] "CMBL"                                    
## [44] "CYP7B1"                                  
## [45] "C8G"                                     
## [46] "IDO2"                                    
## [47] "SLC41A2"                                 
## [48] "DMD"                                     
## [49] "mikado.WCK01-AAH20201022-F8-rsc00002G194"
## [50] "LOC107139914"
FeaturePlot(scSeurat, features = c("PAH","PCK1"))

FeaturePlot(scSeurat, features = c("CYP2E1","FETUB"))

FeaturePlot(scSeurat, features = c("ACACB","ELOVL6"))

FeaturePlot(scSeurat, features = c("POLR2D", "ALB-1"))

Let’s look at some of the key markers that Sonya sent: IL2, IL4, IL6, IL7, IL10, IL12A, TGFB

FeaturePlot(scSeurat, features = c("IL2", "IL4"))

FeaturePlot(scSeurat, features = c("IL6", "IL6-1"))

FeaturePlot(scSeurat, features = c("IL7", "IL10"))

FeaturePlot(scSeurat, features = c("IL12A", "TGFB1"))

It’s pretty hard to see which cytokines are expressed in which samples, so I’m going to try splitting up the samples.

DimPlot(scSeurat, split.by = "tissue_type", label = TRUE) & NoLegend()

Now let’s take a look at the features on this split plot.

FeaturePlot(scSeurat, feature = "IL2", split.by = "tissue_type") # Only in PMAIO cluster 1
## Warning in FeaturePlot(scSeurat, feature = "IL2", split.by = "tissue_type"):
## All cells have the same value (0) of IL2.

## Warning in FeaturePlot(scSeurat, feature = "IL2", split.by = "tissue_type"):
## All cells have the same value (0) of IL2.

FeaturePlot(scSeurat, feature = "IL4", split.by = "tissue_type") # Only in PMAIO cluster 1
## Warning in FeaturePlot(scSeurat, feature = "IL4", split.by = "tissue_type"):
## All cells have the same value (0) of IL4.
## Warning in FeaturePlot(scSeurat, feature = "IL4", split.by = "tissue_type"):
## All cells have the same value (0) of IL4.

FeaturePlot(scSeurat, feature = "IL6", split.by = "tissue_type") # Stronger in endo and mesem PMAIO
## Warning in FeaturePlot(scSeurat, feature = "IL6", split.by = "tissue_type"):
## All cells have the same value (0) of IL6.

FeaturePlot(scSeurat, feature = "IL6-1", split.by = "tissue_type") # No signal
## Warning in FeaturePlot(scSeurat, feature = "IL6-1", split.by = "tissue_type"):
## All cells have the same value (0) of IL6-1.
## Warning in FeaturePlot(scSeurat, feature = "IL6-1", split.by = "tissue_type"):
## All cells have the same value (0) of IL6-1.

FeaturePlot(scSeurat, feature = "IL7", split.by = "tissue_type") # More in PMAIO cluster 4

FeaturePlot(scSeurat, feature = "IL10", split.by = "tissue_type") # More in PMAIO cluster 1
## Warning in FeaturePlot(scSeurat, feature = "IL10", split.by = "tissue_type"):
## All cells have the same value (0) of IL10.

FeaturePlot(scSeurat, feature = "IL12A", split.by = "tissue_type") # No signal
## Warning in FeaturePlot(scSeurat, feature = "IL12A", split.by = "tissue_type"):
## All cells have the same value (0) of IL12A.

FeaturePlot(scSeurat, feature = "TGFB1", split.by = "tissue_type") # Pretty even across all

Cluster 1 looks like a really interesting cluster that has really expanded in the PMAIO experiment. Isolate DE genes for this cluster.

clust1 <- sCVdata_list$res.0.4@DEvsRest$`11`[order(sCVdata_list$res.0.4@DEvsRest$`1`$Wstat,
                                                    decreasing = TRUE), ]
cat(rownames(clust1)[1:50])
## MIB1 FBXL4 PPP3CA OPA1 ST7L CYRIB FILIP1L ARG2 RARB P3H1 TCF4 VPS41 OGA SLC19A2 LOC114104752 AFDN SRP9 ABCC2 TBC1D19 PANK1 FBXO34 ARIH2 OSBPL9-1 SLC22A3 AZGP1 MSH3 SPP2 FBXO25 DLEC1 FPGS TLE4 SLC10A7 B2M HDAC1 TRAM1 FAM135A ANGPTL3 MKKS ABLIM3 TTLL5 B4GALT1 TBC1D12 FKBP8 SCLY PDHX EIF4G2 CNOT10 mikado.WCK01-AAH20201022-F8-rsc00035G807 PACS2 CEP350

Visualise these genes in violin plots.

Interpretation: Now I see that IL2 and IL4 are expressed in more than just cluster 1, but it is nearly always in PMAIO. Cluster 8 is definitely where IL6 is expressed. IL7 is a bit more spread out across treatment types, but it’s definitely present in PMAIO in clusters 4, 5, and 14. Cluster 14 is macrophages, but I’m honestly not sure of 4 and 5 - low signal hepatocytes? IL10 is in clusters 1 and 13 PMAIO. IL12A is a bit more in 6 PMAIO.

VlnPlot(scSeurat, feature = "IL2", split.by = "tissue_type") # Only in PMAIO cluster 1
## The default behaviour of split.by has changed.
## Separate violin plots are now plotted side-by-side.
## To restore the old behaviour of a single split violin,
## set split.plot = TRUE.
##       
## This message will be shown once per session.

VlnPlot(scSeurat, feature = "IL4", split.by = "tissue_type") # Only in PMAIO cluster 1

VlnPlot(scSeurat, feature = "IL6", split.by = "tissue_type") # Stronger in endo and mesem PMAIO

VlnPlot(scSeurat, feature = "IL6-1", split.by = "tissue_type") # No signal

VlnPlot(scSeurat, feature = "IL7", split.by = "tissue_type") # More in PMAIO cluster 4

VlnPlot(scSeurat, feature = "IL10", split.by = "tissue_type") # More in PMAIO cluster 1

VlnPlot(scSeurat, feature = "IL12A", split.by = "tissue_type") # No signal

VlnPlot(scSeurat, feature = "TGFB1", split.by = "tissue_type") # Pretty even across all